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Objectives:  to  model  and  simulate  the  chemical  dynamics  of  0(3P)  reaction  with  hydrocarbon 
monomers,  small  clusters,  liquids,  solids,  and  self-assembled  monolayer  (SAM)  surfaces. 

Accomplishments: 

•  Ab  initio  calculations  on  model  reactions  to  generate  a  database  of  accurate  energies  and 
geometries.1 

•  Least  squares  fit  of  semiempirical  PM3  theory  to  the  ab  initio  database.  The  reparametrized 
potential  is  called  PM3-SRP  (PM3  with  specific  reaction  parameters). 

•  Use  of  PM3-SRP  in  quasiclassical  trajectory  calculations  of  the  dynamics  of  0(3P) 
reactions  with  hydrocarbons. 

Applications: 

o  0(3P)  +  ethane  at  5  eV  collision  energy.2 

o  0(3P)  +  propane  monomer  and  tetramer  at  <  2  eV  collision  energy  (manuscript  in 
preparation).3 

o  Energy  transfer  and  physisorption  of  0(3P)  on  alkane  monolayer  surfaces  (manuscript 
accepted).4 


Introduction.  Reactions  of  the  electronically  ground-state  oxygen  atom  0(3P)  with 
hydrocarbons  are  of  considerable  interest  because  of  their  importance  in  combustion  and 
atmospheric  chemistry,5’6’7  and  in  processing  hydrocarbon  surfaces.8  Extensive  gas-phase  studies 
of  0(3P)  reaction  with  alkanes  (RH)  at  low  collision  energies,  where  the  only  reactive  channel  is 
have  provided  detailed  information  of  product  internal  state  distributions,6’7  scattering  angles  of 

0(3P)  +  RH  OH-+R- 
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the  reaction  products,7  and  rate  constants.9  For  the  0(3P)  +  CH4  and  0(3P)  +  C2H6  reactions 
studied  here,  the  recommended  expressions  for  the  thermal  rate  constants  are,  respectively, 
( 1 . 1 5  x  1  O'15)  F1 56  exp(-4270/T)  cm-*  molecules-'  for  T  =  300-2500  K10  and 

(1.8xl0'31)r6-5exp(-140/T)  cm'3 molecule-1  s']  for  T  =  298-1300  K.11  Some  uncertainties  remain 
in  the  temperature-dependent  rate  constants  as  illustrated  by  the  alternative  expression 
(2.69  x  lO18)!123  exp(-3570/T)  cm  3 molecule-1  sA  suggested12  for  0(3P)  +  CH4. 

0(3P)  is  used  in  chemical  processing  of  hydrocarbon  surfaces  to  make  polymeric 
materials  for  technological  and  industrial  applications.13  Chemical  processing  by  0(3P)  atoms 
may  also  be  detrimental  by  damaging  polymeric  coatings.  This  is  particularly  problematic  for 
spacecraft  in  a  low  Earth  orbit  (LEO).14  Spacecraft  in  LEO  travel  at  a  velocity  of  approximately 
8  km/s,  giving  rise  to  a  relative  translational  energy  of  ~5  eV  for  an  O-atom  striking  the 
spacecraft.  Under  such  harsh  conditions  there  is  considerable  erosion  of  the  spacecraft’s  surface. 
Recently,  several  experimental8-1 5>!  6  and  computational  studies17-18-19  have  been  performed  to 
probe  the  kinetics  and  mechanisms  for  0(3P)  reactions  with  hydrocarbon  surfaces. 

At  high  collision  energies  other  reaction  channels,  in  addition  to  reaction  (1),  are  open  for 
the  0(3P)  atom.  From  ab  initio  calculations,  Massa  and  coworkers20  identified  a  C-C  bond 
rupture  channel,  which  for  ethane  forms  CH3  and  OCH3.  Schatz  and  coworkers 18-1 9  extended 
these  calculations  and  also  identified  a  C-H  bond  rupture  channel;  i.e.  0(3P)  +  RH  — >  H  +  OR. 
For  the  0(3P)  +  C2H6  system  the  threshold  energy  is  approximately  2  eV  for  both  the  C-C  and  C- 
H  bond  rupture  channels.19  For  high  collision  energies,  such  as  5  eV  in  LEO,  it  is  expected  that  a 
large  amount  of  energy  will  be  deposited  in  the  internal  modes  of  the  ethoxy  radicals  by  the  C-C 
and  C-H  bond  rupture  channels,  leading  to  unimolecular  dissociation  of  the  species.  Schatz  and 
coworkers19-21  have  performed  direct  dynamics  simulations  with  the  MSINDO  semiempirical 
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theory  to  study  the  primary  reaction  channels  for  0(3P)  +  C2H6  collisions  at  3.26  eV  to  compare 
with  the  experimental  study21  of  Minton  and  coworkers.  In  this  report  we  describe  a  similar 
direct  dynamics  simulation  at  5  eV  using  a  PM3-SRP  semiempirical  model,  in  which  the 
trajectories  are  integrated  for  a  sufficiently  long  time  to  investigate  the  dissociation  of  the  alkoxy 
radicals  formed  by  the  primary  reaction  channels. 

Multi-reference  ab  initio  calculations }  Electronic  structure  calculations  based  on  multi¬ 
configuration  wave  functions  were  used  to  investigate  the  energetics,  transition  states,  and 
intrinsic  reaction  coordinates  of  a  set  of  prototypical  reactions  for  0(3P)  processing  of 
hydrocarbon  molecules  and  surfaces.  The  specific  reactions  studied  are 


0(3P)  +  ch4 

— - 

OH  +  CH3 

(Rl) 

0(3P)  +  c2h6 

OH  +  C2H5 

(R2) 

OH  +  C2H6 

h2o  +  C2H5 

(R3) 

OH  +  C2H5 

h2o  +  c2h4 

(R4) 

OH  +  C2H5 

h2o  +  ch3ch 

(R5) 

o(3P)  +  c2h6 

ch3o  +  ch3 

(R6) 

0(3P)  +  c2h6 

C2H50  +  H 

(R7) 

For  all  reactions  except  R4,  only  the  lowest  triplet  state  is  examined.  Singlet  spin  states 
are  not  considered.  Reactions  R4-R6  are  possible  secondary  reactions  derived  from  R2.  Reaction 
R5  gives  triplet  ethylene,  3C2H4,  which  is  twisted  90°,  and  R6  gives  triplet  methylcarbene, 
3CH3CH. 

Reactions  R1-R6  were  examined  with  CASSCF,22  second  and  third  order  Rayleigh- 
Schrodinger  perturbation  theory  with  a  CASSCF  reference  (CASPT2,23’24  CASPT323),  and 
internally  contracted  multi-reference  configuration  interaction  (MRCI).25  The  calculations  were 
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performed  using  MOLPRO  2002.326  for  ROHF,  CASSCF,  CASPT2/3,  and  MRCI,  and 
GAMESS  9827  for  the  IRC  calculations  with  CASSCF.  At  CASSCF  and  CASPT2  optimized 
geometries,  single  point  energy  corrections  are  computed  with  MRCI+Q  (MRCI  with  the 
Davidson  correction  for  quadruple  excitations)  and  the  cc-pVTZ  and  cc-pVQZ  basis  sets,  with 
extrapolation  to  the  complete  basis  set  limit  (designated  CBL)  using  two-point  power  law 
extrapolation.28  For  R1-R6,  the  active  spaces  were  (10,10),  (8,8),  (4,4),  (9,9),  (10,10),  (8,8), 
respectively.  Geometry  optimization  and  frequency  calculations  were  carried  out  with 
CASSCF/cc-pVTZ  for  all  reactions,  and  also  with  CASPT2/cc-pVTZ  for  R1-R3. 


Table  1.  CASPT2/CBL  and  MRCI+Q/CBL  values  of  barriers  and  energies  of 
reaction  at  0  K  for  reactions  Rl-R6.a 


Barrier  at  0  K 

best 

estimate 

Energy  of  reaction  at  0  K 

CASPT2 

MRCI+Q 

CASPT2 

MRCI+Q 

exptb 

Rl 

7.9 

10.5 

10  ±  lc 

-0.1 

1.4 

1.6  ±0.2 

R2 

5.3 

9.0 

7  ±  ld 

-4.6 

-1.4 

-5.5  ±  1.4 

R3 

41.3 

46.9 

— 

2.9 

1.1 

-1.9+10 

R4 

0.6 

2.8 

1  -2* 

-18.8 

-17.4 

-21.3  +  1.4 

R5 

1.4 

3.2 

— 

-19.0 

-17.8 

-22.7  ±  4 

R6 

2.0 

3.8 

— 

-10.7 

-10.8 

-14.4  +  5 

aStationary  points  optimized  with  CASPT2/cc-pVTZ  for  reactions  1  and  2, 
CASSCF/cc-pVTZ  for  others. 

bhttp://srdata.nist.gov/cccbdb,  IV.A.l.  Experimental  enthalpies  at  0  K. 
cO.  Roberto-Neto,  F.  B.  C.  Machado,  and  D.  G.  Truhlar,  J.  Chem.  Phys.  Ill,  10046 
(1999). 

^Reference  1. 

^Classical  barrier  from  Y.-Y.  Chuang,  E.  L.  Coitino,  and  D.  G.  Truhlar,  J.  Phys. 

Chem.  A  104,  446  (2000);  ZPE  taken  from  CASSCF(9,9)/cc-pVTZ  frequencies. 


Table  1  summarizes  the  CASPT2/CBL  and  MRCI+Q/CBL  0  K  barriers  and  energies  of 
reaction.  The  best  agreement  (within  experimental  error)  is  found  for  MRCI+Q/CBL  applied  to 
Rl,  the  only  reaction  for  which  the  active  space  preserves  the  full  symmetry  of  reactants, 
transition  state,  and  products,  while  still  including  all  orbitals  that  undergo  significant  bonding 
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changes.  The  MRCI+Q/CBL  0  K  energies  of  reaction  for  R2-R6  are  more  endothermic  than 
experiment  by  3-5  kcal/mol.  A  likely  reason  for  the  discrepancies  is  that,  except  for  R1  and  R3, 
the  CASSCF  active  spaces  do  not  treat  the  C-H  bonds  uniformly  for  reactants,  TS,  and  products. 
Another  source  of  error  in  R4-R6  is  that  MRCI  was  not  able  to  accommodate  reference  spaces 
that  include  all  the  bonds  that  undergo  significant  changes  during  the  reactions. 


Modification  of  PM3  intesrals  and  PM3-SRP  varametrization.  The  potential  V  used  for 
trajectories  is  defined  by 


V  =  v  +v 

r  v  PM  3- SRP  T  r  analytical 


(1) 


where  VPM3_SRP  is  the  PM3  semiempirical  potential  as  implemented  in  the  Mopac  7  package,29 
modified  as  described  below,  and  Vanaiyticai  describes  long  range  O-H  and  O-C  interactions. 

^analytical  is  given  by 


v analytical  =  Z  {[l  +  tanh^.fc  -  r<0))Ja,  exp(-fy;.)+  c/f6] ,  (2) 

/  =  O—C,  O-H 

where  i  refers  to  O-H  or  O-C  interactions,  r(  is  the  distance,  and  Vanajytica]  is  characterized  by  the 
10  parameters  f,  ri0,  at,  bt,  ct  (Appendix,  Table  Al).  VPM3_SRP  is  constructed  with  65  parameters. 
29  of  these  are  conventional  PM3  parameters  for  C,  0,  H  listed  in  Table  A2.  The  remaining  36 
parameters  (Tables  A3,  A4)  are  associated  with  modified  resonance  integrals  (two-center  one- 
electron  integrals)  given  by 

H?=xf(ry)H%ab,  (3) 
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where  Xifb(rij )  is  a  distance-dependent  multiplier  and  is  the  original  PM3  resonance 

integral  between  atoms  i  and  j,  ab  labels  the  overlap  type  (ss,  sp,  pp),  and  ry  is  the  i-j  distance. 
The  multiplier  is  defined  by 


xt(rij)  =  xtmall  +  \{xil^e  -  xts,na„  )£  +  ^(fij  (f*  ~  ^ ))]  >  (4) 


which  switches  from  X^ij, small  at  smaU  values  of  rtj  to  X^ij, large  at  larSe  distance.  O-H 
interactions  were  treated  differently  for  H-abstraction  by  0(3P)  to  form  OH  (OH,  column  in 
Tables  A3,  A4)  and  abstraction  by  OH  to  form  H20  (OH2  column).  This  greatly  improved  the  fit 
to  barriers  and  enthalpies  of  reaction  for  the  two  reaction  types.  Each  trajectory  was  started  with 
OH,  parameters,  and  these  were  switched  to  OH2  parameters  only  in  the  rare  event  that  the 
initially  formed  OH  abstrated  a  second  H-atom.  The  switch  from  OHj  to  OH2  parameters  was 
implemented  when  the  following  distance  criteria  were  met:  O-H  <  1.2  A,  C-H  >  1.5  A,  and  O-C 
>  2.7  A.  Although  the  switch  is  discontinuous,  in  practice  the  tanh  switching  functions  prevented 
sudden  jumps  in  energy  and  gradient  during  the  switch. 

The  65  parameters  of  VPM3_SRP  plus  10  parameters  for  Vanaiyticai  were  determinied  by 
minimizing  the  weighted  sum  of  squared  differences  between  V  and  best  estimates  (ab  initio  or 
experiment)  for  barriers,  energies  of  reaction,  and  geometries  of  stationary  points.  The  fitting 
involves  minimizing  the  sum  of  squares 


»£(«,' -■ E°j  +  y>iZ(Rt- ~ i!,°j  +  ^ »‘A'Z(Ar - A°  j  +  wlX(p, - D°j  +  wl'ZiC,- -  Cl  j 


where  the  difference  terms  involve  PM3-SRP  minus  ab  initio  energies,  bond  lengths,  bond 
angles,  dihedral  angles,  and  Mulliken  atomic  charges,  and  w  are  weighting  factors.  The  Fortran 
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GA  genetic  algorithm  program  of  D.  A.  Carroll  was  used.  The  ab  initio  database  consisted  of  78 
data  points:  26  energies,  20  bond  lengths,  18  bond  angles,  10  dihedral  angles,  and  4  atomic 
charges. 

The  integral  definitions  described  above  have  been  used  in  high  energy  0(3P)  +  ethane 
trajectories2  and  low  energy  trajectories  of  0(3P)  +  propane  monomer  and  tetramer.3 

OOP)  +  ethane  trajectories  at  5  eV  collision  energy.  Full  details  have  been  published,2 
and  a  summary  is  given  here.  Two  parametrization  models  were  examined.  In  Model  1,  the  ab 
initio  database  was  based  on  reactions  R2-R6.  This  is  appropriate  for  low  collision  energies  at 
which  H-abstraction  is  the  only  accessible  reaction.  To  examine  higher  energy  collisions,  single 
reference  ab  initio  calculations  were  performed  at  the  PMP2/cc-pVTZ//UMP2/cc-pVTZ  level  of 
theory.  A  second  parametrization,  Model  2,  was  carried  out  in  which  the  database  included  a 
variety  of  reactions  with  high  barriers  accessible  at  5  eV.  Table  2  gives  a  comparison  of  PM3- 
SRP  Model  2  with  the  ab  initio  data. 

To  calculate  the  classical  trajectories,  Hamilton’s  equations  of  motion  are  integrated  with 
a  combined  4th-order  Runge-Kutta  and  6th-order  Adams-Moulton  predictor-corrector 
algorithm.30  This  was  carried  out  in  VENUS-MOPAC,31  a  combination  of  VENUS32  and 
MOP  AC, 29  which  PM3-SRP  is  incorporated  and  called  whenever  energies  and  derivatives  are 
needed.  There  are  two  criteria  for  terminating  a  trajectory:  (1)  non-reactive  trajectories  are 
terminated  when  the  center-of-mass  separation  between  0(3P)  and  C2H6  is  larger  than  7  A  after 
the  collision’s  inner  turning  point  in  the  0(3P)  and  C2H6  relative  motion;  (2)  reactive  trajectories 
are  integrated  up  to  500  fs  to  monitor  possible  secondary  reactions. 
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Table  2.  PM3-SRP  Model  2  vs  ab  initio  energies  (kcal/mol).a 


Reaction 


PM3-SRP  Model  2 
AEt  AE° 


ab  initio 
AEt  AE° 


Plb 

0(3P)  +  C2H6  OH  +  C2H5 

9.0 

-16.0 

9.0 

-1.4 

P2 

0(3P)  +  C2H6  H  +  C2H50 

34.5 

9.2 

48.3 

11.8 

P3 

0(3P)  +  C2H6  ->  □  d3  +  ch3o 

51.3 

-5.7 

46.9 

1.1 

SI 

OH  +  C2H5  H20  +  3C2H4 

8.9 

-17.3 

4.3 

-17.8 

S2 

OH  +  C2H5  H20  +  3CH3CH 

7.9 

-10.4 

4.9 

-10.2 

S3 

H  +  C2H50  ->  CH2CH20  +  H2 

0.1 

-3.7 

14.5 

1.2 

S4 

H  +  C2H50  ->  CH2CH20  +  H2 

-1.5 

-8.3 

9.6 

-6.0 

S5 

□  □3  +  CH30  — »  □□2QDDd 

6.8 

-12.0 

13.3 

-9.2 

UD1 

c2h5o  ->  ch3  +  ch2o 

16.5 

0.1 

14.4 

5.3 

UD2 

C2H50  ->  H  +  CH3CHO 

9.1 

14.2 

24.7 

6.4 

UD3 

ch2ch2o  ->  ch2  +  ch2o 

15.1 

11.3 

21.8 

12.4 

UD4 

CH3CHO  ->•  CH3  +  CHO 

22.3 

6.5 

9.8 

-3.7 

UD5 

CH30  ->  H  +  CH20 

8.5 

15.0 

20.4 

11.9 

UD6 

CH20  — »  H  +  CHO 

14.2 

20.6 

14.1 

3.4 

aZPE  corrected  (0  K  energies). 

bPl,  P3,  SI,  and  S2  are  calculated  with  MRCI+Q/CBL,  all  others  with 
PMP2/vtz//UMP2/vtz. 


Table  3  gives  cross  sections  for  the  reactive  channels  computed  from  over  50000 
trajectories.  A  substantial  amount  of  the  5  eV  collision  energy  is  deposited  into  the  CH30  and 
C2H50  products  of  channels  3  and  2  and,  as  a  result,  these  primary  products  undergo  secondary 
unimolecular  reactions.  This  leads  to  a  large  number  of  reaction  products  and  small  cross 
sections  for  channels  2  and  3.  Since  channel  1  occurs  by  a  stripping  mechanism,  only  a  small 
fraction  of  the  5  eV  collision  energy  is  deposited  in  the  C2H5  product.  Little  dissociation  of  C2H5 
to  H  +  C2H4  occurs,  as  shown  by  the  large  cross  section  for  channel  1  and  the  much  smaller  cross 
section  for  channel  4.  Because  of  the  importance  of  unimolecular  dissociation  for  the  products  of 
primary  channels  2  and  3,  the  cross  sections  calculated  from  the  trajectories  strongly  depends  on 
the  length  of  time  the  trajectories  are  integrated.  Fewer  product  channels  would  have  been 
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Table  3.  Reaction  cross  sections  of  different  channels. 


Channel 

Products3 

Apparent  reaction 
pathb 

Cross  section,  A2 

Model  1 

Model  2 

1 

OH  +  C2H5 

PI 

7.14 

7.76 

2 

C2H50  +  H 

P2 

0.30 

0.009 

3 

CH30  +  ch3 

P3 

0.025 

0.014 

4 

OH  +  C2H4  +  H 

PI  -»  UD7 

0.082 

0.046 

5 

H20  +  3C2H4 

PI  SI 

0.031 

0.019 

6 

h2o  +  CH3CH 

PI  -»  S2 

0.31 

0.14 

7 

CH3CHO  +  2H 

P2  UD2 

2.06 

1.17 

8 

Ethylene  oxide  +  2H 

P2  ->  UD2  ->  isomc 

0.13 

0.057 

9 

ch2ch2o  +  h2 

P2  S3 

0.021 

<0.001 

10 

3ch3cho  +  h2 

P2  -»  S4 

0.087 

<0.001 

11 

ch2o  +  ch3  +  h 

P2  ->  UD1 

P3  -»  UD5 

0.80 

0.65 

12 

ch2o  +  ch2  +  h2 

P2  -»  S3  -»  UD3 

0.026 

0.042 

13 

CHO  +  CH3  +  H2 

P2  ->  S4  -»  UD4 

0.77 

0.41 

14 

CHO  +  CH4  +  H 

P3  ->  S5  -»  UD6 

0.27 

0.071 

15 

ch2cho  +  h2  +  h 

0.20 

16 

CH2=C=0  +  H2  +  2H 

0.060 

0.072 

17 

C2H20  +  2H2 

0.17 

0.007 

18 

CO  +  ch3  +  h2  +  h 

0.14 

0.93 

19 

CO  +  CH2  +  2H2 

0.084 

20 

CO  +  CH4  +  2H 

0.016 

21 

CH3CHOH  +  H 

P2  ->  isom 

0.061 

22 

3CH2=CHOH  +  h2 

P2  S3  ->  isom 

0.014 

“Ground  state  triplet  potential  energy  surface. 
bNomenclature  for  reaction  paths  is  defined  in  Table  2. 
cIsom  means  the  product  isomerizes. 


observed  if  the  trajectories  were  only  integrated  for  100  ps  instead  of  the  500  ps  calculated  here. 
Indeed,  the  trajectory  cross  sections  may  change  and  additional  products  formed  if  the 
trajectories  were  integrated  for  an  even  longer  time.  Some  of  the  products  observed  at  500  ps 
may  have  sufficient  energy  to  unimolecularly  dissociate  on  a  longer  timescale. 
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OOP)  Reactions  with  Propane  Monomer  and  Tetramer  at  <  2  eV  Collision  Energy 0 
At  low  collision  energy,  H-abstraction  is  the  only  open  channel  for  0(3P)  +  propane,  leading  to 


0(3P)  +  CH3CH2CH3 


OH  +  CH3CH2CH2 
OH  +  CH3CHCH3 


OH  +  CH3CH2CH2 


H20  +  CH3CH2CH 
H20  +  CH3CHCH2- 


OH  +  » -propyl  and  isopropyl  radicals.  Subsequent  reaction  of  OH  +  propyl  on  the  triplet 
potential  energy  surface  can  give  H20  plus  two  products,  triplet  propylene  and  a  triplet  carbene. 
In  our  initial  parametrization  (parameter  set  LI,  where  L  refers  to  low  collision  energy),  the 
energies  in  the  database  included  best  estimate  barriers  and  experimental  0  K  reaction  energies  in 
Table  1.  As  discussed  below,  this  gave  a  hotter  product  OH  rotational  distribution  than 
experiment.33  In  our  current  parameter  set,  L2,  the  database  was  augmented  with  calculations  of 
the  O-H-C  bending  potential  for  0(3P)  +  C2H6  -»  OH  +  C2H5,  using  the  BB1K  density 
functional,34  which  is  optimized  to  reproduce  barriers.  Figure  1  compares  the  O-H-C  bending 
potential  of  BB1K/6-31G*  and  PM3-SRP  with  the  LI  and  L2  parameter  sets.  The  L2  curve  is 
very  close  to  BB1K  at  <  15  kcal/mol  above  the  saddle  point  (minimum  of  the  curve),  and  is 
reasonably  close  at  higher  energies  also.  Table  4  shows  the  PM3-SRP  LI  and  L2  fits  to  reactions 
2-6  used  in  the  parametrization,  as  well  as  the  0(3P)  +  C3Hg  reactions  shown  above.  Barriers  and 
reaction  energies  are  in  reasonable  agreement  with  experiment,  where  available. 

Parameter  sets  LI  and  L2  are  given  in  the  Appendix.  Table  A1  gives  the  parameters  for 
F analytical ,  ecl  2;  Table  A2  lists  conventional  PM3  parameters  accessible  through  the  MOP  AC 
keyword  ‘EXTERNAL’;  Tables  A3  and  A4  give  the  parameters  associated  with  the  distance- 
dependent  integral  multipliers,  eq  3  and  4. 
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Figure  1.  Energy  vs  O-H-C  bending  angle  in  the  transition  state  for  H-abstraction  from 
ethane  by  0(3P),  for  BB1K/6-31G*  and  PM3-SRP  with  parameter  sets  LI  and  L2.  The 
O-H-C  angle  was  varied  at  the  TS  geometry  (O-H-C  =  180°)  without  partial  optimization. 


Table  4.  PM3-SRP  barriers  (AEt)  and  reaction  energies  (AE0),  corrected  for  ZPE, 
compared  to  best  estimates  at  0  K  (kcal/mol). 


Reaction 

□□□□□□ 

param.  set  LI 

□□□□□□ 

param.  set  L2 

Best  estimate 

AEt 

^0 

AEt 

AE0 

AEt 

AE0 

0(3P)  +  C2H6  ->  OH  +  C2H5 

7.4 

-3.7 

7.1 

-5.4 

7  ±  la 

-5.5d 

OH  +  C2H6  H20  +  C2H5 

1.9 

-18.6 

1.9 

-19.7 

1  _2ab 

-21.3d 

OH  +  C2H5  -4  H20  +  C2H4 

4.0 

-15.3 

2.6 

-21.4 

3.2a 

-22.7a 

OH  +  C2H5  ->  H20  +  CH3CH 

1.0 

-12.2 

0.4 

-14.8 

3.8a 

-14.4a 

0(3P)  +  C2H6  ->  CH30  +  ch3 

38.7 

-0.4 

46.9 

14.3 

45. 8C 

-1.9d 

0(3P)  +  C2H6  ->  C2H50  +  H 

44.2 

22.0 

49.7 

41.1 

50. 5C 

8.3d 

0(3P)  +  C3H8  OH  +  h-C3H7 

7.5* 

-4.7* 

6.3* 

1.9* 

— 

-5.4d 

0(3P)  +  C3H8  OH  +  z-C3H7 

5.0* 

-9.5* 

4.4* 

-1.3* 

— 

-8.2d 

*Not  included  in  the  parametrization  database. 

aYan,  T.;  Doubleday,  C.;  Hase,  W.  L.  J.  Chem.  Phys.  120,  9253-9265,  (2004). 
b  Y.-Y.  Chuang,  C.  Coitino,  D.  G.  Truhlar,  J.  Phys.  Chem.  A  104,  446  (2000). 

CD.  Troya,  R.  Pascual,  D.  J.  Garton,  T.  Minton,  and  G.  C.  Schatz,  J.  Phys.  Chem.  A 
107,  7161  (2003). 

dsrdata.nist.gov/cccbdb,  IV. A.  1.,  experimental  reaction  data  at  0  K. 


12 


Extensive  optimizations  of  the  tetramer  yielded  the  global  minimum  in  Figure  2.  This 


structure  maximizes  the  number  of  long  range  attractive  H-H  interactions. 


Figure  2.  Optimized  structure  of  propane  tetramer. 


Approximately  100000  trajectories  were  run  with  the  LI  parameter  set  and  200000  with 
the  L2  set.  Initial  conditions  were  selected  and  trajectories  propagated  using  VENUS-MOPAC,31 
as  described  above  for  0(3P)  +  ethane.  An  important  difference  between  trajectories  run  with  the 
two  parameter  sets  is  the  treatment  of  zero  point  energy.  In  trajectories  run  with  the  LI  set, 
vibrations  are  initialized  with  quasiclassical  normal  mode  excitation.  However,  this  may  not  give 
a  realistic  description  of  zero  point  energy.  Since  normal  modes  are  delocalized  and  are 
initialized  by  superposition  of  randomly  excited  modes,  the  amount  of  energy  placed  in 
individual  C-H  bonds  can  deviate  from  the  desired  value  by  several  kcal/mol,  depending  on  the 
number  of  C-H  bonds  in  the  molecule.  To  produce  a  more  realistic  initial  C-H  vibrational 
distribution,  we  generalized  the  existing  VENUS  code  for  local  mode  excitation  of  a  single 
vibration  to  allow  all  C-H  bonds  to  be  initialized  as  local  modes  instead  of  normal  modes. 
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Quasiclassical  local  mode  excitation  gives  each  C-H  bond  the  desired  amount  of  energy,  and  was 
employed  in  all  trajectories  run  with  the  L2  parameter  set. 

In  all  trajectories,  reactants  were  initialized  with  a  298  K  Boltzmann  distribution  of 
vibrational  and  rotational  levels.  For  propane  monomer,  we  examined  collision  energies  £coll  =  8, 
20,  40  kcal/mol  (0.35,  0.87,  1.74  eV).  For  the  tetramer,  £'coll  =  20,  40  kcal/mol,  and  only 
trajectories  with  zero  impact  parameter  were  examined. 

At  all  collision  energies,  H-abstraction  by  0(3P)  constitutes  98-100%  of  the  reaction,  and 
the  remainder  involves  H20  formation  from  OH  +  C3H7.  Figure  3  shows  opacity  functions  for 
the  monomer  reactions  at  the  three  collision  energies  examined,  for  vibrational  levels  v  =  0,  1,2. 
The  v  =  0  level  is  dominant,  in  agreement  with  experiment.33 


Figure  3.  Opacity  functions  for  reaction  of  0(3P)  with  propane  monomer  at  collision  energies  Ecoh  = 
8.  20.  40  kcal/mol  usina  the  L2  Darameter  set.  Imoact  Darameter  b  is  in  Anastroms. 


A  preference  for  secondary  over  primary  C-H  abstraction  by  0(3P)  is  expected  from  the 
lower  secondary  barrier  (Table  4).  The  results  in  Table  5  are  consistent  with  this.  In  each  case, 
the  percent  secondary  is  above  the  statistical  value  of  25%.  The  preference  for  secondary  may  be 
greater  in  the  tetramer  than  the  monomer,  bu  more  trajectories  would  be  needed  to  confirm  this. 
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Table  5.  Percent  secondary  C-H  abstraction. 


EcoU,  kcal/mol 

Monomer 

Tetramer 

8 

37  ±3 

— 

20 

39  ±3 

43  ±5 

40 

35  ±3 

37  ±4 

Angular  distributions  as  measured  by  solid-angle  differential  cross  sections  vs  cos(vv') 
are  shown  in  Figures  4  (tetramer)  and  5  (monomer).  Here,  v  is  the  reactants’  relative  velocity 
vector  pointing  from  O  to  C3H8,  and  v'  is  the  relative  velocity  vector  of  the  products  pointing 
from  C3H7  to  OH.  As  expected  for  a  nearly  collinear  transition  state,  back-scattering  is  dominant 
at  low  energy.  For  the  monomer,  forward  scattering  increases  with  increasing  collision  energy, 
consistent  with  a  stripping  mechanism. 

Figure  6  compares  the  distribution  of  rotational  levels  j  and  vibrational  levels  v  =  0,  1,2, 
for  trajectories  run  with  parameter  sets  LI  and  L2  at  Eco] ,  =  20  kcal/mol,  and  includes 
experimental  results33  for  v  =  0  (left  plot,  solid  black  circles).  The  important  differences  between 
the  two  sets  of  trajectories  are  that  L2  includes  the  O-H-C  bending  data  in  Figure  1  while  LI 
does  not,  and  that  trajectories  run  with  L2  are  initialized  with  local  mode  C-H  vibrational 
excitation.  The  OH  P(v,j)  distributions  are  clearly  sensitive  to  the  quality  of  the  potential  and  the 
treatment  of  zero  point  energy.  Trajectories  based  on  the  L2  parameter  set,  which  are  initialized 
with  local  mode  C-H  excitation,  give  much  colder  OH  vibrational  and  rotational  distributions,  in 
agreement  with  the  experimental  distribution.33 
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Figure  4.  Normalized  differential  cross  section  (1/ct)[Act/A(cos(vv'))]  for  tetramer  reaction  at  zero 
impact  parameter  and  collision  energy  Ecott  =  20,  40  kcal/mol.  v,  v'  are  initial  and  final  relative  velocity 

vectors.  An  interval  A(cos(vv'))  =  0.2  was  used. 


Figure  5.  Normalized  differential  cross  section  (1/a)[Acr/A(cos(vv'))]  for  monomer  reaction 
at  collision  energy  Ecoi|  =  8,  20,  40  kcal/mol,  where  v,  v'  are  the  initial  and  final  relative 

velocity  vectors.  An  interval  A(cos(vv'))  =  0.2  was  used. 
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EC0||  =  20  kcal/mol 


Figure  6.  Vibration-rotation  distributions  of  OH  derived  from  0(3P)  +  C3H8  monomer  with  Eco!|  =  20 

kcal/mol.  P(v,j)  is  the  probability  of  forming  OH  in  vibrational  state  v  and  rotational  state  j.  Typical  error 
bars  are  shown.  Left:  the  L2  parameter  set  was  used,  and  trajectories  were  initialized  with  local  mode 
C-H  excitation.  Exptl  (black)  from  Ref.  33.  Right:  LI  parameter  set  and  normal  mode  excitation. 


Figure  7.  Vibration-rotation  distributions  of  OH  derived  from  0(3P)  +  C3H8  monomer  with  Ec0]]  =  8 
and  40  kcal/mol  and  the  L2  parameter  set.  Typical  error  bars  are  shown. 
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Figure  7  shows  analogous  P(v,j)  distributions  for  £cojl  =  8  and  40  kcal/mol  using  the  L2 
parameter  set  and  local  mode  C-H  excitation.  The  rotational  distribution  shows  the  expected 
increase  in  the  maximum  of  the  distribution  as  £coll  increases.  Figure  8  gives  the  P(v,j) 
distributions  for  OH  product  from  0(3P)  +  propane  tetramer  with  £coll  =  20  and  40  kcal/mol  with 
the  L2  parameter  set  and  local  mode  C-H  excitation.  Within  the  greater  uncertainty  due  to  fewer 
trajectories,  the  OH  vibrational  and  rotational  distributions  derived  from  0(3P)  +  tetramer  are 
very  similar  to  those  produced  from  0(3P)  +  monomer. 

0(3P)  +  Tetramer  — ►  OH  +  Tetramer(-H) 


Figure  8.  Vibration-rotation  distributions  of  OH  derived  from  0(3P)  +  C3H8  tetramer  with  Eco!|  = 
20  and  40  kcal/mol  and  the  L2  parameter  set.  Typical  error  bars  are  shown. 


Summary  of  OPP)  +  propane.  This  work  presents  a  PM3-SRP  direct  dynamics  trajectory 
calculation  of  0(3P)  +  propane  monomer  and  tetramer  reactions  at  low  collision  energies.  The 
only  important  reaction  at  these  energies  is  H-abstraction  to  form  OH  +  propyl.  Two  parameter 
sets,  LI  and  L2,  were  created  and  extensive  trajectories  were  run  with  each.  In  both 
parametrizations,  the  PM3-SRP  potential  reproduces  the  H-abstraction  barriers  very  well,  and  the 
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reaction  energetics  somewhat  less  well.  Inclusion  of  the  O-H-C  bending  potential  at  the  H- 
abstraction  transition  state  in  the  ab  initio  database  for  the  L2  set  gives  a  PM3-SRP  potential 
energy  surface  on  which  trajectories  reproduce  the  cold  OH  rotational  distribution  observed 
experimentally  (Figure  6).33  Use  of  local  mode  C-H  excitation  gives  propane  the  correct  amount 
of  initial  zero  point  energy,  and  helps  to  reduce  vibrational  excitation  in  the  OH  product.  The 
agreement  with  the  experimental  P(v,j)  distribution  suggests  that  the  methods  employed  here  are 
appropriate  for  reaction  dynamics  involving  general  hydrocarbon  systems. 

Dynamics  of  energy  transfer  and  vhvsisorption  of  OOP)  atoms  in  collisions  with  alkane  surfaces. 

It  is  expected  that  the  reactivity  of  0(3P)  atoms  with  hydrocarbon  materials  will  be 
enhanced  if  a  non-negligible  fraction  of  the  oxygen  atoms  physisorb  on  the  material’s  surface. 
Adsorption  will  provide  more  time  for  0(3P)  to  abstract  an  H-atom  thus  forming  a  hydroxyl 
radical,  or  provide  more  time  for  the  triplet  0(3P)+alkane  system  to  undergo  intersystem¬ 
crossing  to  the  much  lower-energy  singlet  potential  energy  surface  with  the  oxygen  atom 
inserted  into  a  C-H  bond  to  form  R-OH  alcohol. 

Minton  and  co-workers8  have  measured  the  translational  energy  distribution  P(E/)  for 
0(3P)  atoms  scattered  off  300  K  liquid  squalane  at  different  incident  energies  A,  of  5.0,  11.2, 
71.0  and  120.5  kcal/mol.  For  each  of  these  Et ,  there  is  a  300  K  Boltzmann  component  in  P{Ej) 
suggesting  that  a  fraction  of  the  oxygen  atoms  are  physisorbed  and  accommodated  on  (or  in)  the 
liquid  squalane  before  desorbing  with  an  energy  distribution  consistent  with  the  liquid’s  300  K 
temperature. 

We  have  initiated  an  investigation  of  these  experiments  by  first  simulating  inelastic 
collisions  of  0(3P)  with  the  (CioH2iS)n-Au  (s)  1-decanethiol  self-assembled  monolayer 
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(H-SAM).  Due  to  its  homogeneity,  this  is  a  more  straightforward  surface  to  model  than  liquid 
squalane,  and  accordingly  it  is  deemed  to  be  most  productive  to  initially  probe  the  energy 
transfer  and  physisorption  dynamics  for  0(3P)  +  H-SAM  collisions. 


Simulation  Model  and  Potential  Energy  Function 

The  model  used  for  the  alkylthiolate  H-SAM  surface  is  the  same  as  the  one  used  in  our 
previous  studies  of  Ne  atom  scattering  off  this  surface,35’36  with  the  addition  of  periodic 
boundary  conditions.  At  300  K,  the  atoms  of  this  surface  have  non-negligible  thermal 
displacements,  though  the  overall  structure  is  still  well-ordered.  The  following  figure  shows  the 
height  above  the  gold-plane  vs.  time  (0-12  ps)  of  an  H-atom  attached  to  one  of  the  surface’s  top 


O) 

0) 

SZ 


Figure  9.  Plots  of  the  surface’s  H  atom’s  height  above  the  gold 
plane  as  function  of  time  for  ten  different  trajectories. 


methyl  groups.  This  time  dependence  is  plotted  for  ten  different  simulations.  The  plot  shows  that 
the  surface  atoms  sit  about  14.6±0.6  A  above  the  gold  plane. 
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For  the  above  simulations,  the  intermolecular  potential  between  0(3P)  and  H-SAM  is 
considered  critical.  A  QM/MM  model  has  been  developed  for  modeling  chemical  reactions 
between  0(3P)  and  H-SAM,2’35’36*37  and  at  least  one  of  these  models2  may  accurately  describe 
the  intermolecular  van  der  Waals  interaction  between  oxygen  and  H-SAM.  However,  our 
calculations,  which  involve  large  ensembles  of  trajectories  for  many  different  incident  energies 
and  angles,  are  exceedingly  computationally  intensive  and  become  impractical  when  the 
QM/MM  model  is  used  to  provide  the  potential  energy  surface.  Thus,  the  potential  used  here 
consists  of  analytic  functions  for  the  H-SAM  and  an  analytic  function  for  the  oxygen  +  surface 
intermolecular  interactions,  written  explicitly  as  a  sum  of  two-body  interaction  terms  between 
the  oxygen  and  all  the  atoms  of  the  surface.  The  parameters  for  this  two-body  potential  function 
were  determined  by  fitting  three  characteristic  ab  initio  potential  energy  curves,  calculated  at  the 
MP2/cc-pVTZ  level  of  theory,  for  oxygen  approaching  CH4  molecule  in  the  face,  edge  and 
vertex  orientations.  The  potential  energy  curve  for  the  face  orientation  is  shown  in  the  figure 
below. 


O-C  distance  [A] 

Figure  10.  The  MP2/cc-pVTZ  potential  energy  curve  for  0(3P)  +  CH4  in  the  face  configuration.  Red 
points  are  from  the  ab  initio  calculations,  and  the  blue  line  is  a  fit. 
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The  0-H  and  0-C  two-body  intermolecular  potentials  used  to  simultaneously  fit  all 
three  ab  initio  curves  are  of  the  following  form: 

U(r )  =  a  exp(^br)+  c  /  r6 

The  values  of  the  fit  parameters  for  all  three  orientations  are  given  in  the  table  below: 


Table  1 .  Oxygen  atom  +  methane  intermolecular  parameters 


interaction 

a 

b 

c 

[kcal/mol] 

[A-'] 

[mol  kcal1  A-6] 

O-H 

1070.7 

3.191 

-6.56 

O-C 

25106.4 

3.324 

-1044.36 

This  analytic  potential  gives  respectively  an  O  +  CH4  global  potential  energy  minimum  and  O-C 
separation  of  -0.27554  kcal/mol  and  3.499  A  for  the  face  orientation,  of  -0.22756  kcal/mol  and 
3.635  A  for  the  edge  orientation,  and  of  -0.12705  kcal/mol  and  4.086  A  for  the  vertex 
orientation.  A  previous  crossed  beam  experimental  study  of  0(3P)  +  CH4  scattering  gave  an 
average  well-depth  of  -0.22  kcal/mol  and  3.57  A.35  Due  to  a  few  trajectories  in  which  the  oxygen 
atom  penetrates  deep  into  the  surface,  we  have  also  included  O-S  and  O-Au  interactions  which 
are  assumed  to  be  repulsive-only  and  similar  to  O-C  interactions. 

Simulation  Results 


A  Atomic-Level  Dynamics 

The  following  preliminary  results  are  presented  for  simulation  of  0(3P)  +  H-SAM 
collisions  with  a  modest  incident  energy  Et  of  2.30  kcal/mol,  an  incident  angle  <9,  of  60°  to  the 
surface  normal,  and  a  surface  temperature  Tsur/  of  300  K.  The  oxygen  atom  is  assumed  to  begin 
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to  significantly  interact  with  the  H-SAM  surface  when  it  is  at  20  A  above  the  gold  layer  i.e. 
approximately  5.4  A  above  the  terminal  -CH3  groups.  A  similar  criterion  was  used  in  the  study 
of  Ne  +  H-SAM  inelastic  scattering.36  The  minimum  attained  height  of  the  oxygen  atom 
approaching  the  gold  surface  is  found  to  be  always  smaller  than  about  17.5  A.  The  residence 
time  of  oxygen  on  the  surface  is  thus  calculated  as  the  time  period  spent  below  the  20  A  mark.  A 
wide  range  of  residence  times  ris  observed:  The  oxygen  atom  may  bounce  back  as  quickly  as  its 
initial  velocity  allows  it,  with  r<  1  ps,  or  it  may  linger  on  the  surface  for  a  considerable  time,  r> 
20  ps.  Trajectories  also  exhibit  a  wide  range  in  the  number  of  oxygen-atom  hops  on  the  surface. 
Many  trajectories  finish  with  a  single  inner  turning  point,  while  others  have  multiple  turning 
points,  some  in  excess  of  100,  suggesting  a  long-term  physisorption.  The  hops  vary  greatly  both 
in  timing  and  form.  Some  occur  on  the  almost  undetectable  time-scale  (<  0.02  ps)  and  with 
minute  vertical  displacements  (<  0.001  A),  while  others  may  involve  long  times  (>  10  ps)  and 
large  vertical  displacements  (>  4  A).  Depths  of  oxygen  penetration  into  the  surface  may  differ 
greatly.  The  vast  majority  of  oxygen  atoms  sample  only  the  top  few  carbon  layers  of  the  surface, 
but  a  few  penetrate  deep  into  the  surface,  even  hitting  the  bottom  -  the  gold  substrate.  Some  of 
these  latter  oxygen  atoms  appear  to  remain  trapped  inside  the  surface  indefinitely,  i.e.  longer 
than  60  ps  duration  limit  of  the  simulation.  Approximately  2%  of  the  oxygen  atoms  penetrate 
into  H-SAM  below  an  1 1  A  height  from  the  fixed  gold  surface.  We  may  distinguish  five  types  of 
trajectories:  (i)  vast  majority  (98%)  sample  the  top  of  H-SAM  surface  only,  not  penetrating 
deeper  than  13  A,  spending  ~3.5  ps  on  the  surface  and  having  fewer  than  20  hops.  The  other  four 
types  account  for  less  than  1%  each:  (ii)  trajectories  that  remain  on  the  top  of  the  H-SAM  surface 
for  a  long  time  and  experience  numerous  hops,  (iii)  trajectories  that  penetrate  deep  into  the  H- 
SAM  but  do  not  reach  the  Au(s)  bottom  before  scattering  off,  (iv)  trajectories  that  reach  the 
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Au(s)  bottom  but  eventually  resurface  back  and  scatter  off,  and  (v)  trajectories  that  remain 
trapped  near  the  Au(s)  bottom  indefinitely.  The  ~1%  of  very  deep-penetrating  trajectories  are 
associated  with  both  long  lifetimes  and  large  numbers  of  hops  (>20).  Figures  11-16  illustrate 
some  of  the  above-discussed  collision  dynamics. 


Figure  11.  Plot  of  the  oxygen  atom  height  above  Au(s)  for  a  selected  trajectory  with  multiple  hops.  The 
terminal  methyl  groups  are  located  at  -14.6  A.  This  oxygen  atom  remains  near  the  surface  top  at  all 
times. 


Figure  12.  Plot  of  the  oxygen  atom  height  above  Au(s)  for  a  trajectory  with  deep  penetration.  The  atom 
reaches  half-way  down  through  the  surface,  but  manages  to  desorb  within  20  ps. 
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Figure13.  Plot  of  the  oxygen  atom  height  above  Au(s)  for  a  trajectory  which  becomes  caught  in  the 
region  between  Au(s)  and  the  last/bottom  carbon  atom  of  the  alkyl  chain.  The  oxygen  does  not  escape 
during  the  simulation’s  duration.  This  is  an  example  of  very  strong  physisorption. 


Individual  trajectories  are  not  very  informative  nor  representative,  and  thus  it  is  of 
interest  too  look  at  statistics  of  selected  trajectory  characteristics.  In  the  Figures  14-16  are  shown 
probability  distributions  for  some  of  the  significant  gas+surface  collision  aspects  from  a  set  of 
500  trajectories. 
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Figure  14.  Distribution  of  the  minimum  oxygen  height  above  Au(s)  for  the  trajectories. 
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Figure  15.  Distribution  of  the  number  of  hops  that  oxygen  atom  executes  on  the  H-SAM/Au(s)  surface. 
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residence  time  of  O-atom  on  the  surface  [  ] 

Figure  16.  Distribution  of  the  oxygen  residence  time  on  the  H-SAM/Au(s)  surface. 


B  Correlations  Between  Dynamic  Characteristics 

While  dynamic  characteristics  of  individual  trajectories  differ  widely  as  described  by  the 
above  distributions,  there  are  some  important  correlations  between  certain  characteristics  of  the 
0(3P)  +  H-SAM  collisions.  These  characteristics:  the  smallest  height  reached  throughout 
trajectory,  total  number  of  hops,  total  residence  time,  the  final  translational  energy  etc.,  are 
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different  ways  of  describing  already  completed  trajectories,  and  may  have  parallel  trends,  but 
cannot  be  clearly  assigned  any  causal-consequential  relationships  i.e.  dependencies.  It  is  of 
particular  interest  to  determine  how  the  final  translational  energy  Ef  of  the  outgoing  oxygen 
associates  with  various  attributes  of  the  collision  dynamics. 

In  the  following  figures,  some  of  these  correlations  are  presented  in  the  form  of  scatter 
plots.  Figures  17  and  18  show  links  between  Ef  and  the  number  of  hops,  and  between  Ef  and  the 
residence  time,  respectively.  In  collisions  with  no  or  few  hops,  the  oxygen  atom  may  transfer 
anything  from  a  small  to  large  energy  to  the  surface,  while  trajectories  with  many  hops  tend  to 
end  with  departing  oxygen  leaving  behind  the  excess  of  its  initial  energy.  The  Ef  values  for 
collisions  with  Et  =  2.30  kcal/mol  and  with  a  single  inner  turning  point,  vary  from  near  zero 
(elastic)  to  almost  6  kcal/mol,  illustrating  that  in  some  collisions  thermal  energy  of  H-SAM  may 
even  be  deposited  into  the  oxygen  atom.  However,  in  most  trajectories  the  energy  is  of  course 
taken  away  from  the  hotter  oxygen.  Although  the  statistics  is  not  good  enough  with  the  currently 
available  data,  there  are  indications  that  an  increasing  net  energy  transfer  occurs  for  trajectories 
with  more  hops  i.e.  that  oxygen  atoms  that  ramble  more  on  the  surface  are  better  thermalized,  as 
one  would  indeed  expect.  Additional  data  are  nevertheless  needed  to  secure  and  quantify  this 
relationship.  Furthermore,  collisions  with  larger  initial  energy  Et  are  expected  to  show  this 
relationship  more  prominently.  The  trajectories  with  extraneously  large  number  of  hops  are  those 
that  remain  adsorbed  on  the  surface,  and  are  not  included  in  the  figures  below. 
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Figure  17.  Scatter  plot  of  490  trajectories  illustrating  that  the  average  final  oxygen  translational  energy 
(Ef)  (blue  line)  decreases  as  the  number  of  hops  increases. 
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Figure  18.  Scatter  plot  of  the  same  trajectories  as  in  Fig. 17  illustrating  similarly  that  the  final  oxygen 
energy  (Ef)  tends  to  be  small  for  trajectories  with  longer  residence  times. 

For  the  vast  majority  of  trajectories,  there  is  no  noticeable  correlation  between  the  depth 
to  which  the  oxygen  atom  penetrates  the  H-SAM  surface  and  the  oxygen  final  translational 
energy.  It  appears  that  the  energy  with  which  oxygen  is  scattered  off  the  surface  probably 
depends  primarily  on  the  last  few  moments  of  its  residence,  just  before  its  take-off,  and  that  the 
region  where  oxygen  previously  roamed  in  the  surface  i.e.  close  to  the  surface  top  or  deep  inside, 
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is  not  very  significant  except  in  the  sense  that  it  might  tend  to  be  associated  with  different 
residence  times.  There  is  a  weak  positive  link  between  the  maximum  depth  of  penetration 
(minimum  reached  height)  and  the  residence  time.  Indeed,  in  the  few  available  trajectories  with 
deep  penetration,  oxygen  does  not  leave  surface  with  large  Ef .  It  should  be  noted,  however,  that 
all  these  energies  involved  in  present  calculations  are  rather  small.  For  those  few  trajectories 
which  remain  adsorbed  on  the  surface,  there  clearly  can  be  no  Ef  value. 
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Figure  19.  Scatter  plot  of  the  same  trajectories  as  in  Fig. 17  illustrating  little  link  between  the  final  oxygen 
translational  energy  and  the  minimum  height. 

In  Figure  19,  a  correlation  may  be  hidden  or  obscured  by  the  generally  large  scatter  and 
by  the  fact  that  few  trajectories  venture  deep  into  H-SAM.  The  existing  scatter  of  data-points 
represents  the  intrinsic  uncertainty  and  natural  variability  in  trajectory  outcomes,  rather  than  a 
computational  random  error,  of  which  there  is  none. 

Similarly,  there  is  no  noticeable  correlation  between  the  number  of  hops  and  the  greatest 
depth  of  penetration.  In  other  words,  the  oxygen  atom  may  bounce  around  equally  well  in  any 
region  of  the  surface,  in  the  top  as  well  as  in  the  lower  layers  of  the  surface,  as  shown  in  Figure 
20.  Presumably,  the  depth  of  penetration  is  determined  primarily  by  the  relative  surface 
geometry. 
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Figure  20.  Relationship  between  the  number  of  hops  and  the  depth  of  penetration. 


C  Translational  Energy  and  Scattering  Angle  Distributions 

The  distribution  of  the  final  translational  energies  of  the  scattered  oxygen  atom  is  shown 
in  Figure  13,  where  it  is  compared  to  the  300  K  distribution  for  thermal  desorption.  The  initial 
trajectory  conditions  are  the  same  as  above:  Ej  =  2.30  kcal/mol,  =  60°  and  Tsurf-  300  K.  The 
P(Ef)  apparently  has  a  significant  non-thermal  component.  With  an  appropriate  convention,  it 
should  be  possible  to  de-convolute  the  Boltzmann  (BC)  and  the  non-Boltzmann  (NB) 
contributions  when  more  trajectories  become  available  i.e.  when  the  distribution  becomes  better 
resolved. 
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Ef  [kcal/mol] 

Figure  21.  Two  distributions  of  the  final  translational  energies  of  the  scattered  oxygen  atom:  from  the 
trajectory  calculations  (points)  and  statistical  from  thermal  desorption  (line). 


The  probability  distribution  of  the  polar  angle  Of  of  the  outgoing,  scattered  oxygen  atom 
is  shown  in  Figure  14.  The  initial  (incoming)  <9,  and  final  (outgoing)  jangles  are  both  defined  as 
the  angle  between  oxygen  velocity  and  the  surface  normal  long  before  and  long  after  collision, 
respectively.  The  average  outgoing  angle  (Of)  is  47°,  smaller  than  the  incoming  fixed  angle  0{  = 
60°. 
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Figure  22.  Distribution  of  the  scattering  angle  9f  of  the  oxygen  atom. 
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D  Comparison  of  Results  for  Different  Incident  Angles 

The  above  discussion  and  figures  pertain  to  the  atomic-level  dynamics  of  0(3P)  collisions 
at  Ei  =  2.30  kcal/mol  and  6}  =  60°  with  H-SAM.  Additionally,  we  have  performed  analogous  sets 
of  trajectories  with  a  range  of  initial,  incident  polar  angles:  0t  =  15°,  30°,  45°  and  15°.  Each  set 
contains  about  500  trajectories.  Most  of  the  trajectory  properties  have  been  found  to  be  very 
similar,  and  without  any  strong  dependence  on  <9,  except  for  the  minimum  height  of  the  oxygen 
atom  from  Au(s),  suggesting  that  more  slanted  collisions  tend  to  penetrate  slightly  deeper  than 
the  more  orthogonal  ones.  An  interesting  finding  is  that  regardless  of  the  incoming  angle  #  ,  the 
outgoing  angle  Of  is  about  45°,  suggesting  that  the  tilt  and  orientation  of  the  alkane  chains 
determines  the  scattering  angle. 

The  findings  for  all  incident  angles,  and  for  trajectories  that  end  up  with  the  oxygen  atom 
escaping  from  the  surface,  may  be  summarized  as  following:  the  average  final  energy  (Ef  = 
1.5±1  kcal/mol,  the  average  scattering  polar  angle  {Of  =  45°±20°,  the  average  residence  time 
(tres)  =  4±4  ps,  the  average  minimum  depth  of  penetration  (hmif  =  15±2  A,  and  the  average 
number  of  hops  (Nhopf  =  3±4.  In  all  cases,  the  P{Ef  distribution  has  a  significant  Boltzmann 
component,  but  is  still  clearly  non-statistical.  Analysis  is  under  way  to  quantify  the  fraction  and 
characteristics  of  the  non-Boltzmann  contributions  that  can  be  compared  with  the  experiments.8 

The  lack  of  significant  systematic  and  obvious  distinctions  among  trajectory  sets  with 
different  incident  angles  may  be  due  to  several  factors:  these  preliminary  data  may  not  provide 
sufficiently  good  statistics,  or  the  differences  may  be  very  subtle  compared  to  intrinsic  variability 
for  any  simple  comparison.  Also,  the  corrugation  and  complexity  of  the  H-SAM  surface 
structure  may  lead  to  non-trivial  angular  dependence.  Finally,  it  may  be  that  only  the  direct 
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collisions,  with  no  turning  points,  would  preserve  correlation  between  the  incident  angle  and  the 
final  trajectory  characteristics  such  as  the  scattering  translational  energy  Ef .  Further  studies  are 
under  way  to  explore  these  issues,  as  well  as  to  extend  these  simulations  to  higher,  more  relevant 
energies.  It  is  suspected  that  some  of  the  trends,  probability  distributions  and  correlations  might 
become  more  well-defined,  definite  and  distinct  as  one  probes  the  highly  energized  0(3P)  atoms. 

Extensions  to  Of3 P)  +  Liquid  Squalane 

In  simulations  that  are  being  initiated  at  this  time,  the  above  studies  will  be  extended  to 
0(3P)  inelastic  and  reactive  scattering  off  liquid  squalane. 
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APPENDIX 

PM3-SRP  parameter  sets  L2  and  L2  for  0(3P)  +  propane 


Table  Al.  Parameters  for  Vanaiyticai,  eq  2. 


Parameter  set  LI 

Parameter  set  L2 

O-C 

O-H 

O-C 

O-H 

a 

11292.12193 

7808.19694 

11304.51802 

7296.17834 

b 

3.64224 

3.86856 

3.72843 

3.84999 

c 

-509.06505 

-126.28879 

-435.41934 

-109.08288 

f 

16.36333 

20.76417 

19.19875 

27.93941 

r0 

2.52666 

2.96078 

2.49297 

3.01180 

Table  A2.  Conventional  PM3  parameters  of  VPM3_SRP. 


Original  PM3 

Parameter  set  LI 

Parameter  set  L2 

H 

a 

3.356386 

3.553597 

3.633979 

H 

Uss 

-13.073321 

-12.881635 

-12.757258 

H 

Gss 

14.794208 

14.930383 

14.777858 

H 

P-S 

-5.626512 

-6.034970 

-5.950367 

H 

fy 

0.967807 

.958627 

.920356 

C 

a 

2.707807 

2.629564 

2.571816 

C 

Uss 

-47.270320 

-47.073019 

-46.112978 

C 

Upp 

-36.266918 

-35.551186 

-35.647445 

C 

Gss 

11.200708 

12.209925 

11.785205 

C 

Gsp 

10.265027 

10.191357 

10.370604 

C 

Gpp 

10.796292 

10.932203 

10.853004 

C 

Gp2 

9.042566 

8.756955 

8.653818 

C 

lisp 

2.290980 

2.351730 

2.344364 

C 

ps 

-11.910015 

-11.251613 

-11.380953 

C 

Pp 

-9.802755 

-9.413070 

-9.744799 

C 

1.565085 

1.612331 

1.581688 

C 

Cp 

1.842345 

1.818345 

1.797398 

O 

a 

3.217102 

3.196680 

3.191208 

O 

Uss 

-86.993002 

-100.004064 

-97.929636 

O 

Upp 

-71.879580 

-71.349149 

-70.501462 

O 

Gss 

15.755760 

16.406318 

16.447563 

0 

Gsp 

10.621160 

12.384200 

12.052675 

0 

Gpp 

13.654016 

14.292968 

14.129753 

0 

Gp2 

12.406095 

11.281769 

11.051201 

0 

Usp 

0.593883 

.682771 

.669554 

0 

P  s 

-45.202651 

-45.666997 

-46.768261 

0 

P P 

-24.752515 

-24.644967 

-25.440923 

0 

6 

3.796544 

3.949335 

3.818890 

0 

>■ 

2.389402 

2.413209 

2.517966 

V 


Table  A3.  Parameters  of  VPM3_SRP  associated  with  eq  3  and  4,  parameter  set  LI. 


0Ht 

oh2 

CH 

CO 

CC 

HH 

Xs Sij,  small 

0.96737 

0.96737a 

1.00334 

1.03712 

.99433 

.63892 

<ySS..  , 

A  ij,  large 

1.02716 

1.05424 

0.98535 

1.29202 

1.10754 

1.39063 

X^ij, small 

0.95147 

0.95 147“ 

1.00561 

1.01078 

1.05909 

XS^ij,  large 

1.06693 

1.06298 

1.01144 

1.48728 

1.27403 

X^ij,  small 

1.03484 

1.06301 

X^ij,  large 

1.23265 

1.30648 

h 

3.81972 

3.85625 

3.06940 

2.79345 

2.70682 

4.09569 

1.08632 

1.06250 

1.53273 

2.45222 

1.99144 

3.34919 

identical  to  0H[  parameter  (not  independent). 

Table  A4.  Parameters  of  VPM3_SRP  associated  with  eq  3  and  4,  parameter  set  L2. 

OH! 

oh2 

CH 

CO 

CC 

HH 

XSSij, small 

0.944234 

0.944234“ 

1.041469 

1.056987 

1.030936 

0.581784 

XSSij, large 

1.069680 

1.067088 

0.977619 

1.289147 

1.131113 

1.343868 

XSPy,  small 

0.932780 

0.932780“ 

1.006816 

1.017941 

1.039463 

XS^ij,  large 

1.095890 

1.056415 

1.025879 

1.431303 

1.284242 

X^P  ij , small 

1.012786 

1.067242 

X^ij,  large 

1.203813 

1.368468 

fij 

3.703252 

3.811209 

3.081621 

2.783867 

2.725377 

4.095683 

ro.. 

r  V 

1.120038 

1.052181 

1.484811 

2.402244 

1.974958 

3.333669 

a  Identical  to  OHt  parameter  (not  independent). 
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